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1 Introduction 

Computational Fluid Dynamics (CFD) is increasingly being used to both augment and 
create an aerodynamic performance database for aircraft configurations. This aerodynamic 
database contains the response of the aircraft to varying flight conditions and control surface 
deflections. CFD currently provides an accurate and efficient estimate of the static stability 
derivatives, as these involve a steady-state simulation about a fixed geometry. The calculation 
of higher-order dynamic stability derivatives for general configurations and flow conditions 
is more costly however, requiring the simulation of an unsteady flow with moving geometry. 
For this reason the calculation of dynamic stability derivatives using CFD has been limited to 
either estimating the values at a handful of (hopefully) critical points and extrapolating to cover 
the range of interest, or using restrictive approximate methods. The need for more efficient, 
general CFD methods is especially acute as predicting dynamic derivatives with traditional 
methods, such as wind tunnel testing, is expensive and can he error prone. As aircraft designs 
continue to evolve towards highly-maneuverable unmanned systems, high-fidelity aerodynamic 
databases including dynamic derivatives are required to accurately predict performance and 
develop stability and control laws. CFD can provide a key technology for modeling the 
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dynamic performance of these advanced systems, with their extreme rate changes and flight 
conditions. 

The current work presents a novel method for calculating dynamic stability derivatives 
which reduces the computational cost over traditional unsteady CFD approaches by an order 
of magnitude, while still being applicable to arbitrarily complex geometries over a wide range 
of flow regimes. Previous approaches can be broadly categorized as general methods which 
simulate an unsteady motion of the geometry (e.g. a forced osci!lation)[l~5]> or those which 
reduce the problem complexity in some manner with an attendant loss of generality. The 
former methods provide accurate results for arbitrary geometries and flow conditions, however 
they require a time-dependent moving-body flow simulation, which uses roughly an order of 
magnitude greater computational time than a static, steady-state simulation. Even predicting 
one of the roll, pitch, or yaw damping coefficients over the full flight regime is prohibitively 
expensive. Methods which compute the pitch damping using a lunar coning motion[6-8] can 
reduce the unsteady, moving-geometry problem to a static, steady-state computation, albeit in 
a non-inertial reference frame. While these methods provide computational efficiency, they 
are only applicable to longitudinal damping and require approximating (or ignoring) other 
damping coefficients. Weinacht[9] extended this approach to predict pitch and yaw damping, 
however it is only valid for axisymmetric bodies. Weinacht and Sturek[10] demonstrate roll 
damping calculations in a non-inertial frame for finned projectiles, which also reduces the 
problem to a steady-state flow solution, however this approach is only valid at Q: = 0.0°. 
Linearized methods [11, 12] likewise greatly reduce the required computational cost, however 
with a loss of accuracy and a reduced range of applicable flow conditions and/or geometric 
complexity. 

A time-dependent simulation supports a continuum of frequencies up to the limits of the 
spatial and temporal resolution. The primary thesis of this work is that the response to a forced 
motion can often be represented with a small, predictable number of frequency components 
without loss of accuracy. By resolving only those frequencies of interest, the computational 
effort is significantly reduced so that the routine calculation of dynamic derivatives becomes 
practical. Such “reduced-frequency methods” have recently been extended to retain the non- 
linearity of the original governing equations by Hall et al. [13, 14] for application to 2-iJ 
turbomachinery cascades. McMullen et ah [15, 16J followed this work, also focusing on 2-D 
turbomachinery flows. The current implementation uses this same non-linear, frequency- 
domain approach and extends the application to the 3-D Euler equations. 

The current work uses a Cartesian, embedded-boundary method! 17] to automate the gen- 
eration of dynamic stability derivatives. The Cartesian method provides an efficient and 
robust mesh generation capability which can handle an arbitrarily-complex geometry descrip- 
tion. Recently, a method to generate the water-tight surface triangulation required for Cartesian 
mesh generation directly from a CAD representation of the geometry has been developed[18]. 
This, combined with the Cartesian embedded-boundary method provides a robust and auto- 
matic mesh generation infrastructure which can be utilized through the design process. This 
meshing scheme has recently been combined with a parallel, multi-level scheme for solv- 
ing time-dependent, moving-geometry problems, including a generalized rigid-domain motion 
capability [19]. This Arbitrary Langrangian-Eulerian (ALE) rigid-domain motion scheme pro- 
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vides the foundation upon which the current reduced-frequency method is implemented. The 
Cartesian methodology has been demonstrated as an efficient, robust method for automati- 
cally generating static stability derivatives[20, 21], and the current work extends this to include 
the prediction of dynamic derivatives. 

This abstract briefly covers methods for determining dynamic stability derivatives, and 
the reduced-frequency approach, including an 2-D oscillating airfoil example. Demonstration 
examples of calculating the roll damping for the basic finner missile are compared against 
experimental data. The final section discusses topics of future work that will be included in 
the full paper. 

2 Dynamic Derivatives 

The aerodynamic characteristics of an aircraft can be described by the force and moment 
coefficients about the body axes; the axial, normal, and lateral force coefficients (CU, Cyv, Cy), 
and the roll, pitch, and yaw coefficients (Cj, , C n ). In most cases it is sufficient to define 
these coefficients as functions solely of the flight conditions and aircraft configuration, 

Cj = Cj (a,/3, Moo, h, 5i,p, q,r,a, /3^j (1) 

where j represents each of the individual force and moment coefficients, h is the altitude, Si 
represents any configuration-dependent information such as control surface settings, and p y q, 
and r are the rotation rates about the body axes.* Each individual coefficient can be broken 
into two parts: a so-called static portion (subscript s ) which depends only on the non-rotating 
parameters, and a dynamic portion (subscript d) which depends on both the rotational and 
non-rotating parameters. 

Cj = Cj, (a, (3, Moo, h, 5i) + C jd (a, (3, M^, h, 6 it p, q, r, a, /?) (2) 

While the focus of the current work is the calculation of the dynamic portion of the force and 
moment coefficients, a complete evaluation of the static portion is of primary importance. 

In general the aerodynamic coefficients are non-linear functions of all of the independent 
parameters, however in many cases this can be simplified so that a linear superposition of the 
individual effects of each parameter can be assumed, i.e. 

Cj d = C jd (a, /3, Moo, h, 5 u p) + C jd (a, (3, Moo, h,S it q) + ... (3) 

Further, each individual effect can be assumed to be due to a linear variation of that parameter, 
for example the roll variation is given by 

C jd (a, (3, Mooyh, 6i,p) = C jp Ap = (a, {3, M^.h, S t ) A p (4) 


*It is assumed that the rotation rates are suitably non-dimensionalized. 



Notice that these dynamic derivatives arc solely functions of the non- rotating parameters, 
similar to the static coefficients. 

Experimentally, the two tools which are commonly used to provide dynamic data are 
the rotary-balance and forced-oscillation tests. While it is difficult to determine each of the 
individual dynamic derivatives in the general case, as the rotation about the body and wind 
axes are coupled, there is a large legacy of methodology for using data from these tests in 
linearized dynamic models such as described above (cf. Kalviste[22]). The initial focus of the 
current work is to simulate the rotary-balance and forced-oscillation tests using the reduced- 
frequency method, so that the results can be used directly within existing modeling procedures. 
This is seen as a necessary first step; before more complicated uses tor CFD are entertained 
it must become an everyday tool for evaluating dynamic effects in the most common cases, 
similar to the manner it is currently being used to evaluate the static effects. A longer-term 
focus is to develop CFD methods which can compute the dynamic derivatives direcdy in 
the general case, and to extend the methods to provide efficient tools in non-linear flight 
regimes where the traditional methods begin to fail, and obtaining data is extremely difficult 
(cf. Refs. [23, 24]). 


3 Reduced-Frequency Method 


The reduced-frequency method is derived from a general time-dependent scheme. An ALE 
rigid-<iomain motion approach is used in the current application to simulate a forced morion. 
The details of this time-dependent ALE scheme with a Cartesian embedded-boundary method 


rided in [19], and a brief overview is given here. The rime-dependent equations are 


dQ 

dt 


T R (Q) — 0 


(5) 


where t is the physical rime, Q is the vector of conserved variables, and R (Q) is an appro- 
priate numerical quadrature of the flux divergence, ^ j> s f ■ ndS . This work uses an inviscid 
flux vector 

pu n 

f * n = ( pu n u + pn > (6) 

pu n e + pn • n j 


where 


Ur 


(u - u Q ) • n 


is the velocity relative to the moving boundary, and uq is the velocity of the moving domain. 

Following Hall et al. [13, 14], both the conservative variables and R(Q) are assumed to 
be periodic functions of time (with frequency tu), and approximated with a finite Fourier series 


(N- 1)/2 

Q (x, nAt) ^ E 

k—-r-(N—l)/2 


ikunAt 
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(JV— 1)/2 

R (Q, nAt) w Y, *** (Q) e ifeu,nAt 

fc=-(JV-l)/2 

where and R/- are complex Fourier coefficients, and i = ^—l. As a result of this 
approximation Q. can now only support a reduced set of frequencies, namely to and the 
harmonics of u. Since R is a non-linear function of Q, the Fourier coefficient R remains 
a function of Q. Also, since Q and R are real, the Fourier coefficients of the negative 
wavenumbers are complex conjugates of their corresponding positive wavenumber. 

The sampling rate At is chosen so that the functions are periodic over N samples, he. 
co = where T = NAt is the period. The Fourier coefficients can thus be evaluated using 
standard Fast Fourier Transform (FFT) algorithms. Substitution of the Fourier expansions for 
Q and R into the time-dependent equation, Eqn. 5, gives 

ikioQk + Rfc (Q) = 0 (7) 

which form a set of N independent equations due to the orthogonality of the Fourier modes.* 
The solution procedure involves first performing an inverse Fourier transform to construct the 
N samples of Q from These samples are used to construct N samples of R(Q), 

which are then transformed into the Fourier coefficients R^. Equation 7 is then iterated to 
convergence. 

The reduced-frequency approach outlined above has several convenient features. First, it 
can be applied to any set of time-dependent, equations - inviscid, viscous, Reynolds-averaged 
turbulence model, etc. - without requiring any special procedures other than a discrete Fourier 
transform. The same non-linear operator R(Q) from the time-dependent scheme is com- 
puted. Secondly, the same convergence acceleration procedures that are common with static, 
steady-state solvers, such as local timestepping, multigrid, etc., can be utilized to solve Eqn. 7. 
In the current work all of the existing infrastructure from the parallel, multi-level Cartesian 
solver developed for steady-state[25], and unsteady dual-time schemes[19] has been re-used in 
the frequency-domain solver with only minor modifications. 

The cost of the reduced-frequency approach scales as roughly N times the cost of a static, 
steady-state solution, as each iteration requires N evaluations of R(Q). Further, it is required 
to store N copies of each variable in the scheme, which can be prohibitive in 3-D, especially 
using commodity desktop systems. Thus if it requires more than one or two Fourier modes to 
characterize the unsteady behavior, the reduced-frequency approach rapidly loses favor relative 
to solving the time-dependent equations, which require roughly an order of magnitude greater 
effort than a steady-state simulation, but can support a continuum of modes. 

The current approach involves simulating the response to a prescribed periodic motion 
using an invsicid scheme. The simulations of the base state for these flowfields (steady-state 
simulations without a prescribed motion) usually result in a time-invariant flowfield, so that 
all of the unsteadiness in a forced oscillation simulation is due to the prescribed motion. The 
first mode of Q is thus identical to the forcing frequency. Further, in an inviscid simulation no 

*N independent equations for the real and imaginary parts of the positive wavenumbers. The equations for 
the complex conjugate are redundant. 
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physical mechanism exists to transfer energy between modes, so that all of the energy remains 
in the primary mode.* Thus, the response to a forced oscillation computed using just a single 
mode with the reduced-frequency method is equivalent to a full time-dependent simulation for 
an inviscid scheme. This can be seen in Fig. 1, which shows the results of simulating the 
forced oscillation of a transonic NACA 0012 airfoil using three methods: a time-dependent 
simulation, and reduced-frequency simulations retaining one and two modes. After the initial 
transient of the time-dependent simulation, all three simulations are nearly identical, and all 
are in good agreement with the experimental data and capture the hysteresis in the normal 
force variation. This indicates that Q contains just a single mode at the forcing frequency, 
and that including higher harmonics provides no additional information. A comparison of 
numerical timings will be presented for 3-D simulations in the final paper. 



Figure 1: Variation of normal force coefficient with angle of attack for oscillating NACA 0012. (Moo — 0.755, 
a(t) = 0.016 + 2.51 sin (27r62.5t)). The time-dependent simulation includes the initial transient portion of 
the calculation. Experimental data from [26]. 


4 Sample Results 

The reduced-frequency method outlined above is used to calculate the roll damping of 
the basic-tinner missile configuration (Fig. 2) at = 2.5 and angles of attack a ~ 0.0°, 
10.0°, and 20. 0°. The results are compared to a general time-dependent method, and the 
experimental results of Jenke[27j. The time-dependent method computes the rolling moment 

* Numerical dissipation does provide a numerical mechanism to transfer energy, however the numerical 
dissipation is orders of magnitude lower than the physical kinematic viscosity, and hence is not effective. This 
does not imply that a viscous simulation automatically would show a wide energy band, or that these higher 
modes must be resolved to provide an effective estimate of the response to a forced oscillation. 
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on the missile at. constant roll rates p, and then differences the results to determine the 
rolkiamping derivative Ci p . One approach for applying the reduced-frequency method to 
this configuration is to oscillate the body at a fixed amplitude and frequency. Following the 
modeling of dynamic derivatives in Sec. 2, the variation of rolling moment with roll rate is 
assumed to be linear, and hence the roll-damping can be determined from a forced-oscillation 
by a linear regression. The summary of these results are presented in Table 1. The reduced- 
frequency method accurately predicts the roll-damping at all angles of attack at a fraction of 
the cost of the general time-dependent method. 



Figure 2: Basic-finner geometry is a cone-cylinder fuselage with square fins in the + configuration. The cone 
section has a 10° half angle, and the center of mass is located 6.1 diameters from the nose along the longitudinal 
axis of the body. 


a 

Experimental 

Time-Dependent 

Reduced-Frequency 

0.0° 

-17.4 

-17.7 

-17.3 

10.0° 

-19.6 

-19.2 

-18.6 

20.0° 

-21.3 

-21.2 

-19.9 


Table 1: Computed roll-damping C[ p for the basic-finner configuration, M oo = 2.5. 

! 

5 Future Work 

The current abstract contains an overview of the methodology for computing dynamic 
derivatives using a reduced-frequency method. The proposed paper will include a thorough 
validation of the scheme against experimental data for dynamic configurations of interest, 
namely rotary-balance and forced-oscillation test results. These validation cases will include 
i transonic and supersonic results for both complete aircraft and missile configurations. The 

! final paper will also include a summary of the computational cost for the reduced-frequency 

approach compared against the general time-dependent method. 
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